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We propose a realistic scheme to detect the 4D quantum Hall effect using ultracold atoms. 
Based on contemporary technology, motion along a synthetic fourth dimension can be accomplished 
through controlled transitions between internal states of atoms arranged in a 3D optical lattice. 
From a semi-classical analysis, we identify the linear and non-linear quantized current responses of 
our 4D model, relating these to the topology of the Bloch bands. We then propose experimental 
protocols, based on current or center-of-mass-drift measurements, to extract the topological 2"'^ 
Chern number. Our proposal sets the stage for the exploration of novel topological phases in higher 


dimensions. 
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In the past decade, there have been rapid developments 
in the control of coherent quantum systems. These ad¬ 
vances have led to the experimental exploration of fun¬ 
damental quantum mechanical concepts, many of which 
were previously only of theoretical interest. Specifically, 
fascinating artificial systems can now be engineered with 
atoms in optical lattices, where topological order, gauge 
structures, disorder, and interactions can all be tuned 
and probed in novel ways [1-3]. One such important 
recent achievement has been the realization of artihcial 
magnetic helds and topological Bloch bands using time- 
modulated two-dimensional (2D) optical lattices [4-8]. 

Optical lattices can be tailored to realize lattice geome¬ 
tries of varying spatial dimensions dreai = 1,2,3. More¬ 
over, increasing control over atomic internal states [1- 
3, 9, 10] has extended design flexibility as “motion” along 
an auxiliary (synthetic) dimension can be mimicked by 
laser-induced transitions between internal states [11-14]. 
This effectively generates dynamics within synthetic ge¬ 
ometries with (i= (dreai + l) “spatial” dimensions, in ad¬ 
dition to the usual time dimension. This was recently 
demonstrated in experiments, where 2D physics was ex¬ 
plored with atoms in a one-dimensional optical lattice 
with a synthetic dimension [15, 16]. These developments 
naturally open up the possibility of emulating systems 
with higher dimensions d > 3. For example, adding a 
synthetic dimension to a three-dimensional (3D) optical 
lattice would simulate a system in four dimensions (4D). 
Importantly, artificial gauge potentials [2, 3] can natu¬ 
rally be introduced in atomic systems with synthetic di¬ 
mensions [12] in the form of Peierls phase-factors. These 
gauge structures can be finely tuned through the laser¬ 
coupling, allowing experiments to explore topologically- 
nontrivial energy bands in higher spatial dimensions. 

In the fast-expanding field of topological phases of mat¬ 
ter, the energy bands of a system are associated with 
topological indices that are robust to continuous defor¬ 
mations. Nontrivial indices are usually associated with 
interesting boundary phenomena, quantized responses. 


and exotic quasiparticles [17, 18]. In the 2D quantum 
Hall (QH) effect, the quantized Hall conductance is re¬ 
lated to the sum of the 1®* Chern numbers (ICNs) of 
the filled energy bands [19]. In 4D, even more intrigu¬ 
ing quantum Hall phases may exist, as first predicted for 
time-reversal (TR) symmetric systems [20, 21]. In the 4D 
QH effect, energy bands possess an additional nontrivial 
topological index, the 2"'^ Chern number (2CN), leading 
to a nonlinear quantized response [20-22]. Both the 2D 
and 4D QH effects exhibit a variety of exotic strongly- 
correlated phases in the presence of interactions [20, 23]. 
Hence, 4D models with non-vanishing 2CNs have at¬ 
tracted much theoretical attention [21, 24-28]. 

In this Letter, we describe a concrete proposal for re¬ 
alizing the 4D QH effect using ultracold gases in optical 
lattices. We present a semi-classical analysis predicting 
the general transport equations in 4D setups: this in¬ 
cludes a non-linear response, related to the 2CN, but 
also a linear response, associated with an exotic 2D QH 
effect for the current flowing across 2D planes within the 
4D system. We then propose realistic protocols through 
which these responses could be measured experimentally. 

Semi-classical Analysis: We consider particles of 
charge q = —1 [29], moving in a lattice system and 
prepared in a non-degenerate energy band £{k). The 
geometrical properties of the Bloch band’s eigenstates 
ju) are encoded in the Berry curvature 2-form [30], 
n = (l/2)D''‘^(fc)dfc^ A dfc,y, with components = 

t{{dk ,^u\dk„u) — {dk„u\dk,,^u)). The indices yt run over 
spatial coordinates with Einstein summation convention, 
and we set h=l unless otherwise stated. Specifically, we 
are interested in transport under external (weak) electric 
E = and magnetic fields, where B^j,„ = d^Ay—duA^ 
is the magnetic field strength and A = A^e^ is the vector 
potential. We assume that any strong magnetic fields are 
directly incorporated into the Bloch band via the (mag¬ 
netic) eigenstates \u) [31]. Considering for now a sys¬ 
tem of arbitrary dimensions (d > 2), we write the semi- 
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classical equations of motion for a wave-packet as [30-34] 

= ( 1 ) 

k^, = -E^ - (2) 

where r = [resp. k = kfj,e^] denotes the mean po¬ 
sition [resp. quasi-momentum] of the wave-packet. We 
stress that the external fields and are sufficiently 
weak perturbations that, at this level, the wave-packet 
adiabatically follows the band £{k) [35]. Next, we insert 
Eq.(2) into Eq.(l), in order to derive a transport equa¬ 
tion relating the mean velocity to the external helds. To 
second order in the perturbations, this yields [35] 

r^(fc) = + EsB,^^r\k)^^^{k) 

+ (3) 

In addition to the group velocity of the band dk^£, the 
hrst line of (3) includes two important response terms. 
The first of these corresponds to a transverse drift in¬ 
duced by the electric held, which is responsible for the 
QH effect in 2D; the second reveals a non-linear (elec¬ 
tromagnetic) response, which is quadratic in the Berry 
curvature. In this Letter, we are hrst interested in trans¬ 
port of fermions completely hlling the band £{k). In 

this case, neither the group velocity nor the terms in the 
second line of (3) contribute to the response [35]. An¬ 
ticipating that the non-linear term in the hrst line of (3) 
will be crucial when d > 4, as suggested by topological 
held theory [20, 21], we now explicitly consider d = 4. 

4D Transport Equations and Topological Invariants: In 
the semi-classical approach, special care is required when 
evaluating the contribution of all the Bloch states in a 
given Bloch band [30]. An integration performed over 
phase-space involves the usual density of states (l/27r)‘^ 
whenever either the external magnetic held or the Berry 
curvature is absent. However, in the presence of both 
quantities, the physical position and momentum must 
both be distinguished from the canonical coordinates, 
leading to a corresponding modihcation of the density of 
states ]48-50]. In 4D, we hnd that the sum over all Bloch 
states can be accurately replaced by the integral [35] 


3D systems [48]. Combining Eq. (4) with Eq. (3) one 
obtains the transport for a hlled band, i.e. dehning the 
current density yields [35] 

(2^ 1. (5) 

introducing the 2CN [20, 21, 51-53] 

^2 = ^/ DADeZ, 

Jt* 

= ^ f (6) 

47r2 

The transport equation (5) is a central result of this 
work; it was derived from the semi-classical equations of 
motion (l)-(2), without any assumptions about the sys¬ 
tem beyond the non-degeneracy of the Bloch band. The 
only two contributions for the filled band are a linear re¬ 
sponse, reminiscent of the 2D quantum Hall effect [19], 
and a non-linear response proportional to the 2CN. We 
emphasise that the latter effect only takes place in sys¬ 
tems of dimension d > 4, as formally imposed by the 
Levi-Civita symbol. This 4D QH effect was originally de¬ 
rived from a topological-field-theory approach [20, 21] for 
TR-symmetric spinful systems, where the linear response 
vanishes by symmetry [56]. We note that 2CNs were pre¬ 
viously studied from semi-classics in the polarization of 
inhomogenous crystals ]57]. The transport equation (5) 
can also be derived from a Streda-formula approach [35]. 

Realistic Physical System: We now describe our exper¬ 
imental proposal for observing 4D QH physics [Eq. (5)] 
based on cold-atom technology. As is clear from Eq. (6), 
a non-trivial 2CN requires that the Berry curvature 
has non-zero components in “disconnected” planes, e.g. 

ny^ 7 ^ 0. This sets a constraint on the minimal mod¬ 
els displaying non-trivial 4D bands. Inspired by Ref. [28], 
we propose to build a 4D lattice in which the x—z and y—w 
planes are penetrated by uniform magnetic fluxes ^ 1 ^ 2 , 
respectively. This corresponds to two copies of the Hof- 
stadter model [58] dehned in disconnected planes. The 
Hamiltonian is thus of the tight-binding form 




^r+aey'^r 


I i27r$ia:/a t „ , i2TT<i2V/a J 


( 7 ) 

h.c.. 


E 



1 + 


( 4 ) 


-b 


1 




d^k, 


where indices now take values pL = x,y, z, w, V is the sys¬ 
tem’s volume, is the 4D Levi-Civita symbol, and 

where the integration is over the first (magnetic) Bril- 
louin zone (T^). One readily verifies that Eq. (4) gener¬ 
alizes the phase-space integration previously derived for 


where c|, creates a fermion at lattice site r = (x, y, z, w), 
a is the lattice spacing, J is the hopping amplitude 
[Fig. 1(a)]. This model has two main ingredients: (i) x- 
dependent Peierls phase-factors for tunneling along the 
z direction, generating a uniform flux $1 in the x — z 
plane, and (ii) y-dependent Peierls phase-factors for tun¬ 
neling along the w direction, creating a uniform flux 
$2 in the y — w plane. The first ingredient (i) can 
be engineered using the method implemented in Refs. 
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FIG. 1. (a) The 4D lattice in the absence of perturbing fields 
= 0). Atoms move in a 3D optical lattice, with 
hopping amplitude J. A Hofstadter model [58] is realized 
in the x — z plane, with flux 27r$i per unit cell, from x- 
dependent Peierls phase-factors in the hopping along the 2 : 
direction [5-7]. A fourth (synthetic) dimension [11, 12], with 
coordinate w, is accessible at each lattice site, as Raman cou¬ 
pling induces internal-state transitions [35]. A second Hofs¬ 
tadter model, with flux 27r$2, is realized in the y — w plane 
by adding a t/-dependent phase-factor to the internal-state 
transitions, i.e. by aligning the Raman wave-vector along the 
y direction [12, 15, 16]. (b) Energy spectrum £{kj;,ky) for 
$ 1,2 = 1/4, and for many values over the Brillouin zone; 
the 2CN of the lowest (nondegenerate) band is indicated. 


[5-7], i.e., by combining a superlattice (or a Wannier- 
Stark ladder) along the 2 ; direction with a resonant time- 
modulation of the optical-lattice potential [59-61]. The 
second ingredient (ii) requires an extra (synthetic) di¬ 
mension [11, 12, 15, 16], which can be simulated by ad¬ 
justing two-photon Raman coupling to induce transitions 
between different internal atomic states [35]. The Raman 
wave-vector should be aligned along the y direction so as 
to realize y-dependent Peierls phase-factors [Fig. 1(a)], 
as was successfully implemented in Refs. [15, 16]. In gen¬ 
eral, hopping processes along the synthetic dimension are 
inhomogeneous [12, 35], but we verified that this does not 
significantly alter the topological band properties studied 
below [35]. Hence, for the sake of simplicity, we assume 
the 4D system (7) is isotropic with a single hopping rate 
J for all hopping processes. 

The bulk energy spectrum S{k) of the model (7) is 
reminiscent of the two underlying (2D) Hofstadter mod¬ 
els defined in the x — z and y — w planes [28, 35]. Specif¬ 
ically, the spectrum corresponds to the Minkowski sum 
S — E 2 1 .^1 ^ ^xz 7 ^2 ^ ^yw } 7 where Sxz and Syyj de¬ 
note the spectra associated with the 2D Hofstadter mod¬ 
els [35]. This rich 4D spectrum has a lowest band Eiik), 
which is generally non-degenerate and, for suitable fluxes 
of the form — and 4)2 = 1/(72 with ( 71,(72 S Z, can 
be well isolated from higher-energy bands. This happens, 
for example, for fluxes ^12 = 1/4 where the lowest band 
£i{k) has a large flatness ratio A/IF« 3, for bandwidth 
W and lowest bulk gap A«1.3J (Fig. 1(b)). We hence¬ 
forth focus on this choice of 4*1, 2 , as it is relevant to recent 


experiments on the 2D Hofstadter model [5, 7]. 

The bulk spectrum also inherits the topological prop¬ 
erties of the underlying Hofstadter models [28, 35]. 
The lowest band £i{k) is characterized by a non-zero 
2CN [62], which can be factorized as 1^2 = where 

we introduced the ICNs associated with the lowest bands 
of the 2D spectra and respectively [35]. Note 
that the 4D Berry curvature £1^’' and the 2CN V 2 of the 
lowest band can be evaluated numerically, with a 4D gen¬ 
eralization [35] of the method developed in Ref. [63]. 

Response of the System to External Fields: We now in¬ 
vestigate the transport equations [Eq. (5)] with external 
“electric” Ey_ and “magnetic” helds. As the atoms 
are neutral, the “electric” field Ey corresponds here to a 
linear gradient, which can either be magnetically [5, 6] 
or optically [7] created. The perturbing “magnetic” field 
can be generated by engineering additional Peierls 
phase-factors. For the sake of experimental simplicity, we 
choose the weak field By^^ to have a single non-zero com¬ 
ponent, Bzw = —27r$/a^, as is readily achieved by adding 
an extra ^-dependent phase-factor to internal-state tran¬ 
sitions [64]; this occurs in our proposed system if the 
Raman wave-vector is tilted in the y—z plane. To isolate 
the linear and non-linear responses, we align the “elec¬ 
tric” field along the y direction, E = Eyey. The transport 
equations (5) are then explicitly = 0, and 




= E,, 


1 

wr 


n^yd^k = - 


yw 


Ifi 


2'KO?qi 


E,,. 


( 8 ) 

(9) 


Importantly, Eq. ( 8 ) reveals a genuine non-linear 4D-QH 
response along the x direction. In Eq. (9), we have high¬ 
lighted the ICN associated with the y—w plane [65], 
and have used that the flux is 4*1 = 1/(71 in the x—z plane. 
Hence, Eq. (9) shows that an exotic 2D QH effect takes 
place in the y — w plane: defining the Hall conductivity 
in this plane as (Jwy = ~ we find 

awy = -{eVh) i^y^/q,, q,,vy'^ e Z, ( 10 ) 

where the conductivity quantum e^jh is re-introduced. 
Such a “fractional QH effect” is allowed in non¬ 
interacting systems provided the associated 2D plane is 
embedded in a higher-dimensional system [18], as re¬ 
cently shown in 3D topological insulators [ 66 ]. Here, the 
fractionalization of the Hall conductivity reflects the 
change in size of the magnetic Brillouin zone due to the 
flux 4 )i = 1 /( 7 i penetrating the disconnected plane x — z. 

Numerical Calculations: In order to test the predic¬ 
tions in Eqs. ( 8 )-(9), we have numerically calculated the 
current density in a 4D lattice-system filling the low¬ 
est bulk band. The plots in Fig. 2(a) show the quasi¬ 
stationary current density obtained after ramping up 
the field Ey to the final value Ey = —0.2Jla, which was 
chosen so as to limit inter-band transitions [7, 67, 68 ]. 
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FIG. 2. (a) Current density j^{t) after ramping up the “elec¬ 
tric” field to Ey = —0.2J/a. The fluxes were taken 4>i ,2 = 1/4, 
and the Fermi energy Ef = —5J was set in the first bulk 
energy gap. The perturbing magnetic field was chosen as 
—27r4>/a^= -^f5a^. Simulations were performed on a 
8x10 X 10 X 4 lattice, with periodic boundary conditions along 
the X, z, w directions. The component j™ was also calculated 
for $1 = 1/5 (red dotted line). Black lines show predictions 
[Eqs. (8)-(9)] for 1/2 = 1 ^/“= —1. (b) Center-of-mass trajectory 
Xc.m.(^) after ramping up the “electric” field to Ey=Q.2J /a\ 
the non-linear response is in the opposite direction to (a) as 
Ey is reversed. Simulations are for a 12 x 12 x 10 x 4 lattice, 
with periodic boundary conditions along the z, w directions. 
The initial radius of the cloud was ro = 4a in the x — y plane. 
The fluxes are $i _2 = 1/4, with all other parameters the same 
as in (a). The red dotted curve shows the predicted drift, 
u/m. X 16a‘'‘ = — 8 a^U 2 Ey^/'K, for 1/2 = —1- 


Comparing these numerical results to Eqs. (8)-(9) allows 
one to extract the topological indices: we find i/ 2 1.07 
and —1.03, in agreement with the topological band 
structure (i /2 = i/f™ = —1). These simulations also vali¬ 
date the $i-dependence in Eqs. (9)-(10), as we see that 
is indeed reduced by a factor 4/5 when changing the 
flux 4>i = 1/4—7> 1/5 in the x — z plane. 

While current measurements could be performed in 
cold-atom experiments [69, 70], e.g. by generalizing the 
method of Ref. [71], a more accessible probe of the topo¬ 
logical transport would be the observation of the center- 
of-mass (COM) motion of the cloud [7, 68]. The COM 
velocity for this system is = j^Vceii, where is 

given in Eqs. (8)-(9), and V/eii is the (magnetic) unit cell 
volume [7, 68, 72], which is V/eii = 16a^ for ^12 =1/4. 
For a perturbing flux = 1/10, we find that the 2CN 


response should lead to a drift along the x direction with 
^c.m. ~ So/Tb, where Tb = ‘l-n jaEy « 50ms is a typical 
Bloch oscillation time [7]. Such a cloud displacement is 
of the same order as the one reported in the recent ICN- 
measurement [7]. We analyzed the COM displacement 
by simulating a 4D system with open boundary 
conditions in the x,y directions. Numerically, the cloud 
is initially confined in the presence of all fluxes, then the 
confinement is removed and the “electric” field is ramped 
up [68]. From the COM drift in Fig. 2(b), we extract a 
realistic “experimental” value for the 2CN: 1 / 2 0.98. 

Discussion: The COM-response, presented in Fig. 
2(b) for a Fermi system initially filling the lowest energy 
band, is strictly identical to that obtained if the band is 
uniformly, but only partially, filled [35]. This situation 
typically occurs when a thermal gas is loaded into an 
energy band with large flatness ratio [7], such as the low¬ 
est band in Fig. 1(b). Consequently, the 4D QH physics 
discussed here could equally be explored using Fermi or 
Bose gases [7]. We note that other external-field config¬ 
urations could be chosen to fully test the 4D physics in 
Eq. (5). We also stress that an accurate measurement of 
the 2CN requires a large number of internal states [35]; 
this can be solved experimentally by imposing periodic 
boundary conditions along the synthetic dimension, as 
could be realized with an extra Raman coupling [12]. Fi¬ 
nally, we note that the topological gap A « 1.3J is of 
the same order as in [7]; this requires low but achievable 
temperatures to extract the band’s topology (see [7] for 
methods to extract the Chern number in the presence of 
band repopulation, and [73] on heating sources). 

Conclusion: We have proposed a realistic physical 
platform to observe the 4D QH effect in cold-atom trans¬ 
port experiments [7, 8, 67, 68]. We have highlighted the 
minimal requirements for measuring a non-trivial 2CN in 
an atomic system, exploiting Raman transitions between 
internal states [12]. Our proposal sets the stage for the 
future experimental exploration of higher-dimensional 
topological phases. The principles of our set-up could 
also be extended beyond cold atoms, for example, to 
photonic systems [74]. Realizing such 4D systems experi¬ 
mentally will be especially intriguing as they may harbor 
exotic collective excitations [20]. 
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I. CURRENT DENSITY AND TRANSPORT 
EQUATIONS FROM SECOND-ORDER 
SEMI-CLASSICS 

In this Appendix, we show how the 4D current density 
(Eq. 5 in the main text) is derived from the semi-classical 
dynamics of a wave packet. As introduced above, we 
consider particles of charge q = —1 moving in a lattice 
in the presence of an external electric field E = 
and magnetic field strength = d^A„ — d„Afj_, where 
A = Afj^e^ denotes the electromagnetic vector poten¬ 
tial and where indices run over the spatial coordinates 
(Einstein summation convention). We assume that these 
external fields are weak perturbations, which are both 
time-independent and spatially-uniform. We also con¬ 
sider that any strong magnetic fields are incorporated 
directly into the (magnetic) eigenstates |u„) of the un¬ 
derlying Bloch bands with band index n [I, 2]. We stress 
that the following discussion is general and also applies 
to neutral particles, in which case the external fields are 
artificial [3, 4]. 

The wave packet is constructed with a centre of mass at 
position rc = r/e^ and quasi-momentum fe° = A:))e^. Pro¬ 
vided that the external magnetic vector potential A{r) 
varies on a length scale much larger than the spatial ex¬ 
tent of the wave packet, the full quantum Hamiltonian, 
including the perturbing electromagnetic fields, can be 
expanded around Tc as [I, 2, 5-7] 

H + + H" + (I.l) 

where He is the full Hamiltonian evaluated at the centre- 
of-mass position, and H' {H") are first-(second-)order 
gradient corrections in the electromagnetic fields. 

In this semi-classical treatment, we assume that the ex¬ 
ternal fields are sufficiently weak that motion is adiabatic 
and the wave packet remains in a single non-degenerate 
energy band at all times. Eor example, when the electro¬ 
magnetic perturbations are very small, the assumption of 
adiabaticity can be made for a wavepacket constructed 
out of the eigenstates of He [1, 2]. These states corre¬ 
spond to the original (magnetic) Bloch states \un) [up 
to an unimportant phase factor proportional to A{re)] 
[2]. In this case, the resulting semi-classical equations 
capture all effects up to first-order in the external per¬ 
turbing fields. 

Here, we must go a step further as 4D quantum Hall 
physics is nonlinear in the external electric and magnetic 
fields, and so we require a semi-classical treatment that 


is valid up to second order in the perturbations. The 
wave packet should then instead be constructed out of 
the perturbed states [5-7] 

\Un) = \Un) + \u'n) , (1-2) 

where the eigenstates \u'^ are the first-order band-mixing 
corrections from the perturbing term H'. With respect 
to these perturbed states, the semi-classical motion is 
adiabatic even up to second-order in the external fields, 
and we can assume that the wave packet remains at all 
times within the perturbed band n. Hereafter, we drop 
the band index n as well as the subscript c. 

Before proceeding, it is important to note that the ge¬ 
ometrical properties of the perturbed Bloch band are en¬ 
coded in a modified Berry curvature [6] 

- dk^Ak, , 

= (1.3) 

where Ak,, = i{ti\dk,J\u) is the unperturbed Berry con¬ 
nection, containing the effects of any strong magnetic 
fields incorporated into the original Bloch band, while 
A'j.^ = i{u\dk,,\u') -fc.c. is the first-order correction from 
the perturbing electromagnetic fields. We note that the 
correction A'k^ is gauge-invariant, and corresponds to a 
physical shift of the wave packet centre induced by inter¬ 
band mixing from the external fields [6]. 

As derived in Ref. [6], the semi-classical equations in 
dimensions d > 2, valid to second-order in the external 
fields, are 

= (1.4) 

k = -Eu - , (1.5) 

where h = 1. Here, we have also introduced the second- 
order wave packet energy, £{k), which contains the un¬ 
perturbed Bloch band £ (fc) plus corrections from the per¬ 
turbing terms H' H” [6, 7]. 

The above equations of motion can be compared to 
those from usual first-order semi-classics [I, 2]. In this 
latter case, the force equation (1.5) is unchanged, while 
the velocity equation (1.4) is simplified: £l^'^{k) reduces 
to £l>^^{k), and £{k) reduces to £{k) plus the correction 
from H'. As we now show, none of the corrections to 



2 


the band energy nor the Berry curvature affect 4D quan¬ 
tum Hall transport for a filled band; hence, for simplic¬ 
ity of presentation, we have not discussed these in the 
main text. They would, however, have important con¬ 
sequences, for example, in the anomalous Hall physics 
of partially-filled bands, as previously studied in lower 
dimensions [6-8]. 

To proceed, we repeatedly combine (1.4) and (1.5) to 
find the mean velocity: 


d£ 
dk^ 


( 1 . 6 ) 


3f ~ ( 3f 

-h 1 

ak^ i dk^ ' 


d£ 


dku 


d£ 


d£ 


— + + ( ^ + Esn'>^ + , 


dka 


where in the last line, we have neglected terms that 
are above second-order in the perturbing electromagnetic 
fields. We note this approach is similar to that in Ref. [5] 
which derived the contribution from the 2”'^ Chern form 


to the polarisation of an inhomogeneous crystal. 

To calculate the current density of a filled band, we 
sum the mean velocity over all occupied momentum 
states and divide by the volume of the system. In the 
semi-classical limit in 4D, the summation is replaced by 


1 

V 


f d\D{r,k), 

k 


(1.7) 


where V is the volume of the system and where the in¬ 
tegration is performed over the first (magnetic) Brillouin 
zone (T"^). Here, we have introduced D{r,k), the modi¬ 
fied phase-space density of states, which we shall discuss 
in more detail below (see Sec. H). We have generalised 
this function from lower dimensions [9-13], finding in 4D 


D{r, k) 


1 





( 1 . 8 ) 


where all indices now take the values ii = x,y, z, w. 

Combining the above equations and keeping all terms 
up to second order we find the current density 


= f d‘^k[r^D{r,k)] 
Jf* 


Iji 


+ 


d^k 


d£ 

dku 


+ + lEgn^^ + 


f)F 1 f)F \ 


where in the second line we have expanded and grouped 
together terms of a similar form. In this expression, we 
have also replaced the modified Berry curvature H with 
the unperturbed Berry curvature H wherever appropriate 
to ensure that we only keep terms up to second-order in 
the external fields. 

Despite the seeming complexity of (1.10), it is possi¬ 
ble to simplify the current density down to Eq. (5) in 
the main text. We shall now briefly summarise the key 
steps that are required, taking groups of terms in turn, 
beginning with the first two: 




EM^’' 


( 1 . 11 ) 


As can be seen, these terms are very similar to those 
found for a 2D quantum Hall system within usual first- 
order semi-classics [8, 14, 15]. The key differences here 
are the integral over a 4D rather than a 2D Brillouin 
zone, as well as the corrected energy £{k) and the modi¬ 
fied Berry curvature £l{k). However, it can be shown that 


(1.9) 







1 d£ 

2 dkf^ 


B. 


-yi/ 



1 d£ 
64 




( 1 . 10 ) 


both the corrected group velocity dk^E and the Berry cur¬ 
vature correction Di(fc) vanish upon integration over the 
(magnetic) Brillouin zone [6]. The remaining contribu¬ 
tion to the 4D current density is then: 


E, 


/T4 


d'^fc 




( 1 . 12 ) 


which is closely related to the 2D quantum Hall (linear) 
response as discussed in the main text. 

The second set of terms in the current density (1.10) 
are those proportional to both the external electric and 
magnetic fields: 

^4 (1^7 . (1.13) 


This can be simplified by using the antisymmetry of the 
magnetic field strength {B^i, = —B^^) and the Berry 
curvature {W' = —Q£^). By writing out all terms and 
cancelling those which vanish by symmetry, it can be 
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shown that this set of terms reduce to: 

, (1.14) 

where we have introduced the second Chern number for 
a spinless system [16]: 

V2 = -—y [ n A n G Z , 

Jj4 

= (i.is) 

47r2 7 t4 

In the second line, we have written out explicitly the sum 
over the four dimensions for clarity. 

The third set of terms in the current density (1.10) are 
first-order combinations of the external magnetic field, 
the Berry curvature and the group velocity: 


Using only the antisymmetry of the magnetic field 
strength and of the Berry curvature, it can be shown 
that this set of terms sums to zero, giving no contribu¬ 
tion to the integral. We note that this cancellation is by 
symmetry alone, and so holds also for the current density 
of a band with arbitrary filling. 

Combining all these steps, the current density (1.10) of 
a filled band up to second-order in the external fields can 
be simplified to: 

(W L , (L21) 

as stated in the main text. As can be seen, this expression 
does not depend on any of the corrections to the energy 
and Berry curvature that appeared in the full second- 
order semi-classical equations (1.4) and (1.5). 




1 8F 


71 / 




(1.16) 


II. THE MODIFIED PHASE-SPACE DENSITY 
OF STATES D{r,k) IN 4D 


As we have assumed summation over all repeated indices, 
this expression contains many terms; we now group these 
into sets that contain each antisymmetric pair of mag¬ 
netic field strength components. For example, collecting 
all terms with either B^y or Byx{= —B^y), we find: 


B 


d^k 


d£ 


xy 


v (27r)4 I dk^ 


9Fy 


d£ 

Okj. 




d£ 

dk„ 




(1.17) 


where we have used the antisymmetry of the modified 
Berry curvature. As both S{k) and 0(fc) are periodic 
over the Brillouin zone, we can integrate by parts to find: 


f ~ f dn^y dny^ 

Jj4 (27r)4 ^ dki_t dkx ^ dky j ' 

which vanishes due to the Bianchi identity [5]: 

dn^y dnyy^ 

dky dkx dky 


(1.18) 


(1.19) 


The same argument can be repeated for all other pairs 
of magnetic held strength components. As a result, this 
third set of terms makes no overall contribution to the 
current of a hlled band. 

The hnal set of terms in the current density (1.10) are 
higher-order combinations of the external magnetic held, 
the Berry curvature and the group velocity: 


/T4 


d^k 


BF 1 f)F \ 

akn 2 ak-y 


1 d£ 
64 dka 




. ( 1 . 20 ) 


As discussed above and in the main text, to calcu¬ 
late the transport equation for a hlled band in the semi- 
classical limit [Eq. (1.9)], we integrate the mean veloc¬ 
ity over the Brillouin zone, taking into account a modi- 
hed phase-space density of states [9-12]. We now discuss 
the interpretation of this factor and derive its extension 
[Eq. (1.8)] to four dimensions. 

In the absence of the Berry curvature and/or the ex¬ 
ternal B held, the phase-space density of states is a con¬ 
stant factor, given by l/(27r)'^. Classically, the indepen¬ 
dence of the density of states from the phase-space co¬ 
ordinates is a consequence of Liouville’s theorem, which 
states that the phase-space volume element is conserved 
under time evolution [9]. However, when both the mag¬ 
netic held and Berry curvature are non-zero, it can be 
shown that Liouville’s theorem applies to the volume el¬ 
ement AU = ARAK of the canonical position R and 
momentum K, but not necessarily to the volume element 
AU' = ArAk of the physical position r and momentum 
k [9-11]. This is because the relationship between the 
canonical and physical variables is, in general, strongly 
affected by the presence of the magnetic vector potential 
and the Berry connection. 

In the simplest case, when we have a magnetic held 
but no Berry curvature, it is well-known that the phys¬ 
ical momentum is modihed by the Peierls substitution 
k = K — A{r), while the physical position remains equal 
to the canonical position [17]. If instead, we have Berry 
curvature but no magnetic held, the physical position be¬ 
comes r = R + A{k), while the physical and canonical 
momenta remain the same [8, 18, 19]. This is because 
the Berry curvature acts as an artihcial magnetic-held in 
momentum space [19, 20]. Consequently, in the presence 
of both Berry curvature and a magnetic held, neither the 
physical position nor momentum remain equal to their 
canonical counterparts. Instead, we must apply gener¬ 
alised Peierls substitutions for both variables which con- 
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tain many additional nontrivial terms; this is discussed 
in detail in Refs. [6, 8, 11, 12]. 

Returning to Eq. (1.7), we see that the violation of Li- 
ouville’s theorem for the physical variables means that, 
in general, the phase-space density of states with respect 
to r and k is not a constant. Hence, when we replace the 
sum over physical momenta by an integral, we include 
the factor D(r,k) which can depend on both the phys¬ 
ical position and momentum. We emphasise again that 
while the violation of Liouville’s theorem may at first 
seem surprising, the dynamics are of course still Hamil¬ 
tonian because of the preservation of Liouville’s theorem 
for the canonical variables [10, 11]. 

To understand D{r, k) further, there are several equiv¬ 
alent approaches that have been discussed in the litera¬ 
ture [9-12, 21]. Firstly, as in Ref. [9], we could consider 
the time evolution of the physical volume element given 
by: 


1 dAV' _ d d ■ 

AW dt ~ 


(H.l) 


Substituting in from Eq. (1.4) and Eq. (1.5), this can be 
solved analytically for AV'{r,k). The modified density 
of states is then chosen to ensure that D{r, k)AV'{r, k) 
remains constant at all times, i.e., that the change of 
the volume element under time evolution is compensated 
for by the change in the density of states [9]. In this 


sense, it is natural to call D{r,k) the modified phase- 
space density of states. 

A second interpretation of D{r,k) is as the Jacobian 
of the transformation from the canonical to physical vari¬ 
ables multiplied by a factor l/(27r)'^ [10-12, 21]. In this 
picture, we see that the integral of the velocity in Eq. (1.9) 
is constructed first for the canonical variables, where Li¬ 
ouville’s theorem applies and the phase-space density of 
states is l/(27r)‘’*. Then the integral is transformed by 
the Jacobian into the physical variables. In the absence 
of either a magnetic field or the Berry curvature, it can 
be shown that the Jacobian is equal to unity and the 
transformation is unimportant [10, 13]. 

This Jacobian description also leads to simple deriva¬ 
tion of D{r, k) from classical Hamiltonian dynamics [10]. 
As this is the most straightforward calculation for 
D{r,k), we apply this method to find the modified 4D 
phase-space density of states for Eqs. (1.4) and (1.5). 
In classical Hamiltonian dynamics, the Euler-Lagrange 
equations can be written as : 

^ . ( 11 . 2 ) 

where h is the classical Hamiltonian, ^ denotes the phase- 
space physical coordinates r and k collectively, and uja/s 
is a symplectic matrix. Re-writing Eqs. (1.4) and (1.5) in 
this form, we have in full: 


/ 

0 

Bxy 

-Bxz 

Bxw 

-1 

0 

0 

0 

\ 


f 

X 

\ 


f 

Ex 

Ey 

Ez 

Ew 

de 
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Bxy 

0 

-Byz 

Byw 

0 

-1 

0 

0 




y 






Bxz 

Byz 

0 

Bzw 

0 

0 

-1 

0 




Z 






Bxw 

Byw 

Bzw 

0 

0 

0 

0 

-1 




w 






1 

0 

0 

0 

0 

n^y 


^xw 




hx 




dkz: 



0 

1 

0 

0 


O 

ny^ 

^yw 




ky 




d£ 

dky 



0 

0 

1 

0 


-ny^ 

0 





kz 




di 


V 

0 

0 

0 

1 

_^xw 

_^yw 


0 

) 


[ 

kyj 

J 


[ 

dkz 

d£ 

dkyj 

) 


where we identify the classical Hamiltonian a.s h = S + 
E^r^. Then, it can be shown that the Jacobian of the 
transformation from the canonical to physical variables 
can be simply calculated as the square-root of the deter¬ 
minant of the symplectic matrix [10]. In 4D, this is: 

D{r, k) = \/det(w) (H.4) 



as stated above and in the main text. In less than four 
dimensions, the last term vanishes due to the Levi-Civita 


symbols, and what remains is 

^2D/3D(r, fe) = 1 + , (11.5) 

which is the previously known result for two and three 
dimensions [9, 10]. 

Finally, we note that while the above discussion was 
framed in terms of classical variables, these issues can 
also be approached from a quantum perspective, for ex¬ 
ample by re-quantizing the semi-classical equations [9] 
or by beginning from the full quantum Hamiltonian and 
performing a semi-classical diagonalization [12, 21]. 
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III. THE STREDA-WIDOM FORMULA 


The Streda-Widom formula is a simple thermodynamic 
relation to calculate the quantum Hall response of an 
insulating material [22-24]. While this was originally 
derived for the linear 2D quantum Hall effect, it has 
since been discussed in the context of nonlinear physics in 
Ref. [25], where the magnetoelectric coupling of a three- 
dimensional topological insulator was studied. We now 
show that a direct application of the Streda-Widom for¬ 
mula appears to also capture the physics of 4D quantum 
Hall transport. 

The bulk of a system is insulating provided that the 
chemical potential lies in a band gap. Under the appli¬ 
cation of an electric field, local equilibrium is established 
with the electric field Ei, imposing a gradient of the chem¬ 
ical potential /x across the system [24]. If the electric 
field is sufficiently weak, the bulk remains insulating as 
the chemical potential stays in the gap. However, the 
boundaries of the system can support surface currents, 
which are a function of the magnetic moment evaluated 
at the local chemical potential [23, 24]. From thermody¬ 
namics, it can be shown that the current density of the 
surface currents is given by [24]: 


J = 


dn 




T,B 


E„ 


(HI.l) 




where M'^‘^ is the magnetic moment, n is the particle 
density and T is the temperature. In 2D, it has been 
shown that this current density may be quantised in units 
of the first Chern number; this is the standard quantum 
Hall effect [9, 22-24]. 

We now apply this formula in four dimensions; then 
the density of particles is 


As can be seen from (HI.2), the only dependence of 
the particle density on a perturbing magnetic field, at 
a constant chemical potential, comes from the modified 
phase-space density of states D{r,k) [9]. Then consider¬ 
ing, for example, the surface current density along the x 
direction under an applied electric field Ey, it is straight¬ 
forward to show that: 


J = 


dn 


dB, 


= E„ 


yx 


It* 


(27r 


Ey 

T,y 


F, B 

^y-^zw 


V2 

i2ny 


(HI.3) 


where we have used that the corrections Hi vanish upon 
integration over the whole Brillouin zone, and have only 
retained terms up to second-order in the fields. We have 
also applied the definition of the second Chern number 
for a spinless system (1.15). This procedure can be re¬ 
peated for all other configurations, recapturing the gen¬ 
eral result for the 4D current density stated above (1-21); 
this suggests the central role of the modified 4D phase- 
space density of states in the physics of the 4D quantum 
Hall effect. 


IV. SPECTRAL PROPERTIES AND CHERN 
NUMBER DECOMPOSITION 

In this Appendix, we elaborate on the spectral details 
of the model studied in the main text [cf. Eq. (7)]. We 
are considering a 4D model that is a direct sum of two 
copies of the 2D Hofstadter model [26]. For emphasis, let 
us rewrite Eq. (7) in another form 

H = -J (IV.l) 

'w,x,y,z 


n= f d^kD{r,k)<d{£p — e), (III-2) 

Jt* 


as extended from lower dimensions [9]. Here, Q{£f — s) 
is the distribution function of fermions, which uniformly 
fill the band up to Fermi energy ep- We assume that the 
Fermi energy lies in a band-gap, such that the band is 
completely filled. 


where 


LHofstadter/ r /f,\ 

^a./3 —^a,l3,j+l,S^a,l3,j,S 


-A 


(IV.2) 


f27r<J>7 t 


^a,/3,7,(5-|-l^a,/3,7,(5 h.C. , 


and here we have assumed that the lattice spacing a = 1. 

The spectrum of each Hofstadter model is obtained by 
solving the 2D Harper equations 


^XZ{kx7 kz )4’x{k X-) kz) 

^yw {ky, kyj)\y(ky, ) 



E'^^(j)x+i{kx,kz) + e "'^^(jjx-iikxjkz) - 2^,j;(fca;, cos(27r$ia: + kz)j, 

e ^ \ypi {ky , kw) ^Xy — 1 i^y 5 ) 2xy (A:^, k^o^ cos(27r4?2y -i- k^ ) 'i , 


(IV.3) 

(IV.4) 


where x = {!,...,gi}, y = {I,..., 92 } for $1 


Pi/qi and $2 = P 2 / 92 , and where Pi,2,9i,2 are inte- 
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(a) 



kx kx ky ky, 


FIG. IV.1. (a) The energy spectrum £{kx,ky), as shown in 
the main text, of the 4D model for <I>i ,2 = 1/4, and for many 
values over the Brillouin zone. The 2CN of the lowest 
band is indicated. This 4D spectrum is a Minkowski sum of 
the spectra S^zik) and Sywik) of the 2D Hofstadter models 
in the x — z [resp. y — w] plane, shown in (b) on the LHS 
[resp. RHS]. The first Chern numbers of the lowest band for 
each 2D Hofstadter model is indicated. 


gers; here the eigenfunctions satisfy the periodicity con¬ 
ditions 4>x{kx,kz) = 4>x+qiikx,kz) and Xy{ky,kw) = 

Xy+q 2 iky,ky,). For example, the spectrum of each 
2D Hofstadter model with flux = 1/4 is seen in 
Fig. IV. 1(b). 

The energy spectrum of the full 4D model, as seen in 
the main text and here in Fig. IV. 1(a), is a Minkowski 
sum of its two Hofstadter constituents, i.e., 

£{k) = {Ei + E2 \Ei e £xz{kx,kz),E2 G £yw{ky,ky,)}. 

As a result, filling up the lowest band, up to a chemical 
potential p, within the first gap of £{k), corresponds to 
a Minkowski sum over the energies of the lowest band of 
each Hofstadter constituent. Thus, when evaluating the 
2'^'^ Chern number, 1/2 [cf. Eq. (6) in the main text], the 
integral can be split such that each 2D Berry curvature 
is integrated up to the lowest band of the corresponding 
Hofstadter spectrum, and we obtain the result appearing 
in the main text 1^2 = x , where we introduced the 
1®* Chern numbers 

= n^^ik^,k.)dk^dk,, (IV.5) 

XTT Jj2 

associated with the lowest bands of the 2D (Hofstadter) 
spectra £xz and £yw, where the integral is over the 2D 
(magnetic) Brillouin zone (T^) for momenta ky^,k^. 


V. NUMERICAL EVALUATION OF THE 
BERRY CURVATURE AND 2ND CHERN 
NUMBER 

In this Appendix, we briefly describe a numerical 
method allowing for the computation of the Berry curva¬ 
ture £l>^‘'{k) in 4D lattice-systems, i.e. yL.v = x,y,z,w. 
In particular, this efficient method offers a powerful 
method to numerically evaluate the 2'^'^ Chern number 
of 4D Bloch bands, defined in Eq. (1.15) above. We 
emphasize that the algorithm described here is a direct 
generalization of the method developed by Fukui et al. 
[27] for 2D lattice systems; we thus refer the reader to 
this reference for technical details. 

We start by considering a general 4D lattice sys¬ 
tem. The momentum-space representation of the lattice 
Hamiltonian is written as H{k), where k — {ky} is the 
quasi-momentum, taking values within the first (mag¬ 
netic) Brillouin zone. 

For the sake of simplicity, we will consider that 
the spectrum £{k) associated with H{k) displays a 
non-degenerate, well isolated, lowest energy band [e.g. 
Fig. IV. I (a)]. Generalization of the method to the case 
of non-degenerate bands directly follows from the deter¬ 
minant formula [Eq. (16)] in Ref. [27]. 

The numerical method relies on a discretization of k- 
space: H{k) is solved on a fc-space lattice, with lattice 
sites denoted fc^. The ground eigenstate \u{kL)) is nu¬ 
merically calculated at each lattice site, which allows us 
to build a link variable 

D'^(fcL) = (u(fci)lu(fei + e^))/J\f, (V.I) 

where A/” is a normalization constant. In lattice gauge 
theory [27], such link variables are directly related to the 
curvature through a plaquette product. Here, the (con¬ 
tinuous) Berry curvature £l{k) = {k)dkyAdki, is 

approximated by [27] 

n^^''{kL)=i \og\u^{kLW{kL + e^^) (V.2) 

X [U^{kL+e'')]-^[U’^{kL)]-^]. 

The 2'^'^ Chern number is eventually calculated using the 
discretized expression 

^2 = ^ (V.3) 

We applied this simple algorithm to our model 
[cf. Eq. (7) in the main text], for the flux configuration 
chosen in the main text. The Hamiltonian matrix H{k) 
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is associated with the (Harper-like) eigenvalue equation 
£{k)ux,y{k) = -jl^e'^'"^Ux+i,y{k) + Ux-i,y{k) 

+ (V.4) 

- 2ux,y{k) [cos(27r<I)ia; -I- kz) + cos(27r$22/ + ku,)] \ , 


In the following, we consider a single band, filled uni¬ 
formly with a filling factor p; the case p = 1 corresponds 
to that of an exactly filled Bloch band. In this partially- 
filled-band configuration, the current density in Eq. (5) 
is simply modified by the constant filling factor, 

(VL2) 


where x = {l,...,( 7 i}, y = {l,...,g 2 } for $i = Pi/qi 
and <I >2 = P 2 lq 2 , and where pi, 2 ,qi ,2 are integers; here 
again the eigenfunctions satisfy periodicity conditions 
'^x.yik^ — '^x+qi,yik^ and Ux^y{k^ — ^a:,y+q2(^)’ 

For $ 1^2 = ±1/95 the fc-space lattice was discretized as 

27r 

ky. = ■^ith nx,y = 1 ,..., 91; nz,w = 1 ,..., gOl . 

Stable and accurate results were obtained for 91 > 3, 
namely, for fc-space lattices of size larger than 3 x 3 x 
3q X 3g. We first verified that the only non-zero compo¬ 
nents were £l^^{kx, kz) and which 

is in agreement with the configuration of the fluxes 
and d> 2 , which are penetrating the x—z and y—w planes, 
respectively. We also verified that the lowest energy band 
of the model was indeed associated with a non-zero 2'^'^ 
Chern number given by ^2 = [cf. Sec. IV]. In par¬ 

ticular, the numerical method yields a remarkably precise 
value V 2 = vf^ = = —1. for <I>i 2 = 1/4, see Fig. IV.1, 

and91>3[27]. 


VI. PARTIALLY-FILLED BANDS AND THE 
CENTER-OF-MASS DISPLACEMENT 


which directly follows from a generalization of the calcu¬ 
lations in Section I of the Supplemental Material. 

We relate the current density to the centre-of-mass 
(COM) velocity of a cloud of atoms by noting that the 
former is given by jp = t>tot/fo) while the latter is 
W.m. = Vtot/Ntot, for a total particle number TVtot and a 
total velocity Vtot- Combining these definitions, we see 
that Vc.m. = jp/n where n is the particle density. For the 
configuration of lattice model and perturbing fields con¬ 
sidered in the main text, the density of a partially filled 
band is simply n = p/14eii- Then the COM velocity 
is 


= 


Nt, 


p/Vcell 


= j^VceU- 


(VI.3) 


Hence, unlike the current density, the COM velocity is 
independent of the filling factor p. The center-of-mass 
measurement described in the main text could therefore 
be equally carried out with systems that only partially 
fill the lowest (flat) band [e.g. Fig. IV.1(a)] as with those 
that completely fill the band. This indicates that the 
4D QH physics we have described will be experimentally 
accessible with both Fermi and Bose gases. 


In this Appendix, we briefly discuss transport for a 
partially-filled energy band. In this case, the current den¬ 
sity in Eq. (1.9) is replaced by 

jji= f d^kp{k)f^D{r,k), (VI.l) 

where p(fc) is a filling function, which can be determined 
by band-mapping techniques [28]. As in the main text, 
we still restrict ourselves to the situation where only the 
lowest band is occupied; generalization to multi-band 
configurations is straightforward, see Ref. [28]. 

In general, p(fc) may be an arbitrary function. How¬ 
ever, in the experimentally-relevant case where the oc¬ 
cupied band has a large flatness ratio A/FF [e.g. as in 
Fig. IV.1(a)], both the Berry curvature £l^''{k) « 
and the filling function p(fc)« p are uniformly distributed 
over the magnetic Brillouin zone. The latter is because 
the temperature is generally large compared to the band¬ 
width, yet still can be made small (or of the order) of the 
band gap, W <C A:bT ^ A. In fact, the interplay between 
temperature, inter-particle collisions and the periodic- 
driving of the lattice can even lead to a uniform pop¬ 
ulation of dispersive bands, as was recently observed in 
Bose gases loaded into Hofstadter bands [28]. 


VII. OPEN BOUNDARY CONDITIONS ALONG 
THE SYNTHETIC DIMENSION 

In this Appendix, we analyse the possibility of measur¬ 
ing the 2"*^ Chern number without engineering periodic 
boundary conditions along the synthetic dimension (w). 

Indeed, the numerical results shown in Fig. 2 [main 
text] were obtained in a configuration where V = 4 in¬ 
ternal states were used to form a closed synthetic di¬ 
mension; from the experimental point-of-view, this ne¬ 
cessitates the further coupling of those states delimiting 
the edges in the w direction (which we loosely denote as 
r(; = l and w = N), see Ref. [29]. Here, we instead con¬ 
sider a finite 4D lattice with open boundary conditions 
along all directions, including the synthetic one. We note 
that current experimental setups using ^^^Yb atoms [30] 
can potentially generate a synthetic dimension with up 
to V = 6 internal states (from the / = 5/2 nuclear spin 
manifold). To go further still, the ground-state mani¬ 
fold of atoms could be used; this has a total angular 
momentum F = 9/2, potentially allowing experiments to 
work with N^IO internal states. 

The Chern-number measurement introduced in the 
main text relies on the detection of the center-of-mass 



(COM) displacement. This method measures a bulk 
property, and as such, it remains accurate as long as the 
edges of the system do not affect the motion of the cloud. 
In particular, the method breaks down if the cloud sud¬ 
denly hits the sharp edges at w = 1 and w = N. There¬ 
fore, even if the system is initially well-localized around 
the central site (w) « A^/2, corresponding to preparing 
the gas in a small number of internal states, this method 
would require the use of a large number of internal states 
fV»1 to be applicable in practice. 

We now investigate this finite-size effect more quanti¬ 
tatively. Let us consider the relevant setup configuration 
illustrated in Fig. 1 [main text], for which the COM 
velocity is given by the expressions 

v^„,=fxl6a‘^ = 8a^Ey^/TT, (VII.l) 

<m. = r X lQa^ = 2a^Ey/Tr, = 0. 

Here, we used Eqs. (8)-(9)-(VI.3), and have taken into 
account the fact that 1^2 = this configura¬ 

tion as the fluxes are $ 1,2 = 1/4. In the following, we 
also set the perturbing flux to $ = 1/10, and the elec¬ 
tric field strength Ey = 0.2J/a, as used in the main text 
[see Fig. 2]. From the COM velocity (VII.l), one can 
readily determine the time t* at which the cloud touches 
the edge in the w direction. Let us suppose that the 
cloud is initially prepared around the central site (with 
w-coordinate r“ « aN/2), with a spread of about 1-2 lat¬ 
tice sites in the w direction; this corresponds to initially 
populating 1-2 internal states. Then the critical time is 

« « (y -1) « Y (f - 1) (vn.2) 

where we have reintroduced H. This indicates that the 
duration of the COM-drift measurement should not ex¬ 
ceed the critical time t*^10h/J when working with TV « 4 
internal states. 

Next, we remember that the 2”'^ Chern number is 
extracted from the COM displacement along the x 
direction [Fig. 2(b)]. It is thus relevant to estimate 
the typical displacement Xc.m.^t*) performed by the 
cloud during this critical time. Using Eq. (VII.l), 
we find Xc.m.i^*) ~ ®/2 for N « 4. Hence, in the 
absence of engineered periodic boundary conditions, the 
2CN-measurement would necessarily require single-site 
resolution imaging; we note that a differential measure¬ 
ment would show a COM displacement of about one 
lattice site, 2xc.m.(C) « a. In order to overcome this 
technical difficulty, one should either use more internal 
states (we find 2a;c.m.(C) who for VwlO, which could be 
realized by populating the ground-state manifold of 
atoms), or more preferably, engineer periodic bound¬ 
ary conditions [29], as already presented in the main text. 


VIII. SYNTHETIC DIMENSIONS AND 
ANISOTROPY EFFECTS 

In this Section, we investigate the effects of anisotropy, 
which are inherent to Raman-induced synthetic dimen¬ 
sions for cold atoms [29]. To address this question, 
we consider the relevant configuration of Ref. [29], in 
which Zeeman sublevels [mi?) in the hyperfine ground- 
state manifold (with total angular momentum E) are 
coupled by two Raman lasers, with frequency difference 
UJ 1 —UJ 2 = wr and wave-vector difference ki — k 2 = fcR. 
The energies of the Zeeman sublevels uimp are shifted 
by a (real) constant magnetic field, and we assume that 
the resonance condition a;R= Wm_F-i-ii® fulfilled for 
all mp (i.e. we neglect detuning and quadratic Zeeman 
effects; see Ref. [29]). This atom-light-coupling configu¬ 
ration is described by an effective Hamiltonian [4] 

i/eff = , F± = Y ± iPy, 

where we introduced the ladder operators acting on the 
ground-state manifold 

F+Iitif) = gF,mp\mF + 1), (VH.3) 

PF{F + 1) — mpimF + 1), (VH.4) 

and where Hr is the generalized (two-photon) Rabi fre¬ 
quency [4]. This effective Hamiltonian can be inter¬ 
preted as a tight-binding Hamiltonian describing tunnel¬ 
ing along a synthetic dimension spanned by the 2F +1 
internal states, 

2 F +1 

i/eff = ^ ct+iC„ + h.c., (VH.5) 

W — 1 

Jw(w) = Hr^ F{F -I- 1) - + 1), (VH.6) 

where creates a particle at lattice site w, is the tun¬ 
neling amplitude, and where exp (T/sr • x) is the Peierls 
phase-factor exploited in the main text. In Eq. (VH.6), 
we introduced the coefficient TOu, =tof corresponding to 
the internal state \mF) living at site w] in the following, 
we use the convention according to which the first site 
ty = I hosts the lowest-energy state with mp = —F, the 
next site w = 2 hosts the state with mF = —F + 1, and so 
on. 

Atomic-lattice systems based on synthetic dimensions 
are potentially anisotropic [29], due to the fact that the 
tunneling rates Jw in Eq. (VH.6) explicitly depend on the 
w coordinate, through the coefficients gF,Tap- A method 
by which anisotropy can be limited in the full lattice (in¬ 
cluding real and synthetic dimensions) consists in tuning 
the Rabi frequency to the value 

2F 

^K = J/{g), (5) = X] 9F,7npl2F, 

m—1 

where J is the hopping amplitude along the real di¬ 
mensions (we assume that tunneling rates associated 
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with real dimensions can be made homogeneous and 
isotropic). In this case, the tunneling amplitudes are 
given by 


Jyj{w) = J F{F + 1) - m^{m^ + l)/{g)^ « J Vw. 

(VIL7) 

We now investigate the tunneling rates using the con¬ 
figuration defined in Eq. (VII.7). For F = I, we have 
9F,mF = for all Zeeman states, so that the system 
is genuinely isotropic in this case: = J- In the case 

F = 3/2, we obtain the following tunneling rates Jw{w) 
between the four lattice sites: 


, 0.95J „ I.IOJ „ 0.95J . 

re = I- w = 2 - w = 3 - w = 4. 


In the case F = 9/2 (ten sites), we find 


T 0.71J „ 0.95J „ 1.09J . 

w = 1 -^ w = 2 - w = 3 - w = 4 

. 1.16J 1.19J „ 1.16J „ 

w = 4 - w = 0 - w = 0 - w = I 


1.09J^ 

w = 7 - w = 


0.95J „ 0.71J 

-s> w = 9 -s> w = 10 


The effects of anisotropic tunneling rates Jw was al¬ 
ready studied in detail in Ref. [29], for an atomic lattice 
made of one real and one synthetic dimension; this sys¬ 
tem exhibited the quantum Hall effect and was thus char¬ 
acterized by a band structure with nonzero Chern num¬ 
bers. In that study, it was shown that the anisotropy due 
to Jw only slightly modified the bulk bands, and in par¬ 
ticular, that this imperfection preserved all topological 
gaps (even in the case of large F); see the Supplemen¬ 
tary Material of Ref. [29]. 

We have performed a similar study for the 4D system 
considered in this work, and we verified that the topo¬ 
logical band structure remains remarkably robust against 
anisotropy. In the ideal case of an isotropic 4D lattice 
[see Eq. (7) in the main text], the bulk energy spectrum 
depicted in Fig. 1(b) displays a non-degenerate lowest 
band. This lowest bulk band, which is associated with 
a non-zero 2nd Chern number, is separated from higher- 
energy bands by a gap of size A = 1.32J and is char¬ 
acterized by a reasonably small bandwidth W = 0.43J. 


These values for the gap and bandwidth were obtained 
for an isotropic 4D system of size 40 x 40 x 40 x 40 with 
periodic boundary conditions. In order to study the ef¬ 
fects of anisotropy, we have computed the bulk energy 
spectrum for the anisotropic 4D system, including the 
realistic tunneling rates described by Eq. (VII.7), for 
the cases F = 3/2 (four sites) and F = 9/2 (10 sites); 
the corresponding system sizes are 40 x 40 x 40 x 4 and 
40x40x40x 10, respectively. Similarly to the 2D case 
discussed in Ref. [29], we obtain reasonably small devia¬ 
tions from the ideal (isotropic) configuration. In the case 
F = 3/2, we obtain a topological gap of size A = 1.39J 
and a bandwidth W = 0.36J. These values are to be 
compared with A = 1.43J and W = 0.31,/, which cor¬ 
respond to the values obtained for an isotropic system 
with four lattice sites along w. This negligible modifica¬ 
tion of the topological band structure indicates that the 
simulations presented in the main text [Fig. 2], which 
were realized for an isotropic system with four lattice 
sites along w, satisfactorily capture the physics of a re¬ 
alistic (anisotropic) lattice configuration. For F = 9/2, 
we obtain a gap A = 0.90J [and W = 0.84J], which is 
still reasonably large as far as the Chern-number mea¬ 
surement is concerned. From this analysis, we conclude 
that all topological properties described in the main text 
for the ideal isotropic lattice are well preserved in the re¬ 
alistic anisotropic lattice. We emphasize that the (2nd) 
Chern-number measurement entirely relies on the topo¬ 
logical properties of the populated band and on the size 
of the bulk gap A, which are both well preserved under 
anisotropy effects. 

Finally, we point out that periodic boundary condi¬ 
tions can be applied along the synthetic dimension, us¬ 
ing an additional Raman transition [29]; this is relevant 
when dealing with a small number of internal states; see 
Section VII. Importantly, the tunneling rate associated 
with the closing link, i.e. Jw(w = 2F + 1), can be con¬ 
trolled independently, which can be exploited to further 
reduce anisotropy. The results discussed above in this 
Section were obtained by setting Jw(w = 2F+1) = J. 
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